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D . We compute shear viscosity and flavor diffusion coefficients for ultra-relativ- 

Q \ istic gauge theory with many fermionic species, N{^>1, to leading order in 

1/iVf. The calculation is performed both at leading order in the effective 
coupling strength g 2 Nf, using the Hard Thermal Loop (HTL) approximation, 
£T) • and completely to all orders in g 2 N{. This constitutes a nontrivial test of how 

well the HTL approximation works. We find that in this context, the HTL 
\ approximation works well wherever the renormalization point sensitivity of 

the leading order HTL result is small. 

o 

I. INTRODUCTION 

Many problems in early universe cosmology, and in heavy ion physics, require thermal 
field theory. Perturbation theory was developed for thermal field theory decades ago [1]. 
However, the application of perturbative techniques to problems sensitive to long time scales, 
or to low frequency or wave number plasma excitations, runs into subtleties which have made 
it very difficult to compute interesting time dependent phenomena even at leading order in 
the gauge coupling. 

The physics of Fourier modes with frequency and wave number (u, k) ~ gT (with g the 
gauge coupling and T the temperature), the so called "soft" degrees of freedom, cannot be 
treated even at leading order in g without the resummation of a class of diagrams, called 
the Hard Thermal Loops (HTL's) [2]. This problem was elucidated by Braaten, Pisarski, 
Frenkel, Taylor, and Wong, who showed that, once the HTL diagrams have been re-summed 
into the propagators and vertices of soft excitations, one may make a loopwise expansion 
which is now an expansion in g. This breakthrough has made possible the calculation at 
leading order of a number of time dependent plasma properties, such as particle damping 
rates [3] and hard particle energy loss [4]. It has also allowed the development of an effective 
theory for even more infrared fields, relevant to baryon number violation in the standard 
model [5]. 

Another problem of interest in thermal field theory is to compute very long time scale 
properties of the plasma. Since the most infrared properties of the plasma admit a hydro- 
dynamic description, this requires, besides understanding of the thermodynamics, the cal- 
culation of transport coefficients (shear and bulk viscosity, electric conductivity, fermionic 
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number diffusion). The transport coefficients can be formally related to zero frequency and 
momentum limits of correlation functions of physical observables, and require for their com- 
putation, at leading order in the coupling, the resummation of an infinite class of diagrams, 
which for scalar field theory has been shown to be ladder graphs [6]. 

To date, the only leading order calculation of a transport coefficient in relativistic field 
theory is for one component scalar field theory [6]. Despite a substantial literature [7-19], it 
is only recently that transport coefficients have been computed correctly in a gauge theory 
to leading order in the logarithm of the coupling [20]. This is partly because there are 
new complications in hot gauge theories; unlike scalar field theory, hydrodynamic transport 
coefficients are sensitive at leading order to soft physics 1 . Hence their calculation requires 
control of both long time scales and soft momenta. 

Clearly it would be valuable to compute transport coefficients to full leading order in the 
coupling within a gauge theory, even if only in some simplifying regime which rendered the 
calculation easier. Such a calculation will require resummation of hard thermal loops. It 
would be even more valuable to compute transport coefficients beyond leading order in the 
coupling. Aside from the intrinsic interest of the calculation, such a calculation would also 
permit a test of the quality of the HTL approximation. The best measure of the usefulness 
of a perturbative expansion is to see how fast the series converges as one goes to higher and 
higher order. So far, in applications where HTL resummation is necessary, even the best 
results are at leading order in the coupling 2 . Therefore it is hard to know what confidence 
to assign to a result computed at leading order using hard thermal loops. 

The purpose of this paper is to compute transport coefficients, at leading order and to all 
orders in the coupling, in a toy theory where it is possible to perform the calculation. The 
toy theory is QED or SU(iV c ) QCD, with a large number N f of massless fermions, N { ^> N c , 
iVf > 1. We treat g 2 as very small, and will work only to leading nontrivial order in this 
quantity; but the t'Hooft-like coupling g 2 N can either be expanded in perturbatively, or 
treated as 0(1). As we will see, at leading order in 1/iVf, nonabelian effects are unimportant, 
so we will generally use language and normalizations of QED. 

QED with Nf ^> 1 is still a very nontrivial theory, with substantial similarities to realistic 
QED or QCD. Unlike a scalar theory, the theory is derivative coupled. HTL corrections are 
needed for soft gauge boson lines, and represent a rich set of plasma physics phenomena, 
such as Debye screening and Landau damping. At weak coupling, scattering cross sections 
show the Coulombic divergence, cut off by screening effects; there is also interesting collinear 
physics with an analog in full QED or QCD (see Sec. VB). Yet the theory is enough simpler 
that we can treat it to all orders in g 2 N { . Naturally, this is because it is missing some of the 
physics of QED or QCD. For instance, bremsstrahlung is absent at leading order in 1/Nf, 
even at large g 2 Nf. The added complications in QED and QCD may mean that the HTL 
expansion is worse behaved in those theories, so our results on the convergence should be 
considered as optimistic estimates for the series convergence in realistic theories. 



1 For the case of bulk viscosity, even in the scalar theory there is sensitivity to soft physics [6]. Bulk viscosity 
is more complicated for other reasons, and will not be considered further here. 

thermodynamic quantities are a notable exception, where for instance the pressure is known to fifth order 
in the coupling [21]. 
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FIG. 1. All bubble graphs, left, and the resulting self-energies when they are cut, right. 

The large N { theory has a peculiar structure; while there are an enormous number, 
O(Nf), of fermionic degrees of freedom, there are only 0(1) gauge field degrees of freedom 
(for QCD we treat N c ~ 1). Since the vertices of the theory always involve the gauge boson, 
the rate at which a fermion interacts is actually very small; there are only 0(1) degrees of 
freedom for a fermion to couple to, and g 2 ~ l/N { is very small. The gauge bosons, on the 
other hand, couple to all of the different fermions. Though the interactions are individually 
weak, g 2 N{ is not small, so for the gauge bosons the interactions are significant. 

This structure leads to key simplifications which make the calculation of transport co- 
efficients tractable without an expansion in g 2 Nf small. Since the fermions only couple to 
0(1) gauge boson degrees of freedom, the fermionic self-energy is 0(g 2 ) ~ 0(1/N{) <C 1. 
This behavior is clear at one loop, but it is also true to all orders in loops. The easiest way 
to see this is to draw all bubble diagrams which appear in the theory at O(iV f ), which are 
illustrated in Fig. 1. Since any self-energy can be generated by cutting one line in a bubble 
diagram, the fermionic self-energy to the order of interest can be obtained from this set of 
bubble diagrams. 

The smallness of the fermionic self-energy means that fermions are almost free particles. 
In particular, the imaginary part of the fermionic self-energy is ~ 1/iVf for all momenta. 
Therefore, there is a sharp quasi-particle mass shell for fermions, even if g 2 N{ is paramet- 
rically 0(1). Fermions propagate freely over very long distances, undergoing occasional 
collisions with mean separation ~ Nf/T. For this reason a kinetic theory description of the 
fermions is reliable, with 0(g 2 ~ 1/iVf) (and therefore negligible) corrections. 

A kinetic description of fermions only, without including gauge bosons as kinetic degrees 
of freedom, is also sufficient. This is so because fermions dominate the observables of interest 
in this paper, the stress-energy T^ u and all fermionic number currents, simply because there 
are O(iVf) times as many fermionic as gauge degrees of freedom. 

The other key simplification is that, at leading order, the fermion and gauge boson self- 
energies each have a very simple structure. The gauge boson propagator is, at 0{Nf) but 
to all orders in g 2 iVf, given by re-summing a one loop fermionic self-energy, 

D'l = - 11^(1 loop) . (1.1) 

This should be clear from Fig. 1. The fermionic self-energy is gotten by putting the re- 
summed gauge boson propagator into the self-energy diagram shown in Fig. 1. At leading 
order in a 1/iVf expansion the fermionic lines needed to compute IL^ in Eq. (1.1), and the 
fermionic line appearing in the fermionic self-energy, are free theory lines. 

These simplifications will be sufficient so that we will be able to solve for the transport 
coefficients within kinetic theory, without uncontrolled approximations. Our method for 
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solving for the transport coefficients follows very closely our previous work [20], except that 
we will not need to make a leading log approximation of scattering processes, but will be 
able to treat them in a complete way (at leading order in 1/iVf). Some integrals need to 
be done numerically and the momentum dependence of the departure from equilibrium has 
to be modeled with a several parameter Ansatz, but for both approximations the error in 
the treatment can be made arbitrarily small with sufficient numerical effort; modest effort 
drives relative errors to 0.1%. 

The large Nf theory has a technical problem, which is that it does not exist, because the 
gauge boson propagator exhibits a Landau pole at finite Q 2 . However, provided (^Landau T 
(say, at least 40T), this problem is only manifested at energy scales which the thermal 
bath probes exponentially rarely. It is therefore irrelevant to the thermal physics we study. 
Avoiding Landau pole physics restricts somewhat the range of g 2 N we can consider, but 
it is not a severe restriction; the perturbative expansion is still not guaranteed to be well 
behaved at the largest values of g 2 N { which satisfy the above criterion. 

An outline of the paper is as follows. We discuss how the transport coefficients are to 
be computed via kinetic theory in Sec. II. The most complicated feature of the kinetic 
theory is the collision integral, which is discussed in Sec. Ill, with some details relegated 
to two appendices. Our exact (at leading order in l/Nf) results for transport coefficients 
appear in Sec. IV. Then, Sec. V discusses how the HTL approximation can be used in the 
calculation to obtain the result to leading order in g 2 N { , using a simpler scattering matrix 
element. It presents a comparison of the leading order result with the exact result. Finally, 
there is a conclusion, Sec. VI. A very brief summary of the conclusion is that, at relatively 
weak coupling, the leading order calculation works extremely well; but at couplings relevant 
for even the highest energy heavy ion collisions currently being considered, a full leading 
order treatment is at best a factor of 2 estimate, because of severe renormalization scale 
uncertainty. The renormalization point suggested by dimensional reduction turns out to do 
surprisingly well at large Nf] but it is difficult to extrapolate, whether this will occur more 
generally. 

II. TRANSPORT COEFFICIENTS AND KINETIC THEORY 
A. Transport coefficients 

The equilibrium state of a plasma, or any system, is determined by the densities of 
all conserved quantities. For the theory we consider, massless QED or QCD with many 
Dirac fermion species, the conserved quantities are the energy-momentum tensor, fermionic 
particle number, and the charges under SU(A c )xSU(A c ) flavor rotations. In the grand 
canonical ensemble these are determined by the timelike inverse temperature 4-vector (3^ 
and various chemical potentials // (with s a combined species and spin label). Here (3^ 
determines both the temperature 1/T = (3 = ^J— (3^(3^ and the velocity of flow Ui = (3i/ (3q 
of the fluid, and the fi s determine the values of the conserved particle numbers. [Note that 
we use a ( — h++) index metric throughout, so timelike 4-vectors have negative squares.] 

The hydrodynamic regime is the regime where the densities of conserved quantities vary 
slowly enough in space that the local state at each point can be described up to small correc- 
tions by the equilibrium configuration, but with spacetime dependent (3^ and fi s . Transport 
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coefficients are defined in terms of the departure from equilibrium, due to the gradients of 
conserved quantities, and are related to entropy production and the damping away of the 
spatial inhomogeneity. In particular, viscosities are defined by the difference between the 
equilibrium and actual stress energy tensor. If in equilibrium and in the (local) rest frame 
u(x) = the stress tensor is given by T^- = SijV (V the pressure, a function of /3 M and the 
fi s ), then when u has gradients, T^- is given to lowest order in the gradients of u by 

r 2 _ i 

Tij = S i:j V - 77 [Viiij + VjUi - -<Jij-Vi«iJ - (SijViUi , (2.1) 

with i] the shear viscosity and ( the bulk viscosity. (Both are functions of (5 and the /i s .) 
Evaluating bulk viscosity is quite complicated because it involves the physics of particle 
number changing processes [6]; we will not treat it here, but estimate parametrically that 
( ~ NfT 3 , times a dimensionless function of g 2 Nf, in this theory. 

The other transport coefficient we consider is number diffusion. The diffusion constant 
for a conserved current is defined by the constitutive relation (in the frame u = 0) 

ji = -DVdo. (2.2) 

If we choose to consider a U(l) theory where all fermions have the same U(l) gauge charge, 
the total fermionic particle number is constrained, in equilibrium, to be zero. However this 
restriction does not apply for nonabelian gauge theories, or for the other conserved particle 
numbers. There are therefore a large number of fermionic number densities with diffusion 
constants we can consider. In general we should permit the diffusion "constant" to be a 
matrix in the flavor space of conserved currents. However, to leading order in Nf all the 
diffusion constants prove to be the same, so we will ignore this complication. In the U(l) 
theory, the diffusion of the total fermionic number density is related to the behavior of 
infrared electric fields and turns out to give the electrical conductivity; this is discussed 
more in [20]. 



B. Kinetic theory 

The stress tensor arises predominantly from fermionic degrees of freedom, simply because 
there are O(Nf) of them, and only 0(1) gauge degrees of freedom. Naturally, fermionic cur- 
rents and densities are determined entirely by the fermionic degrees of freedom. Furthermore, 
because (? 2 ~ 1 /A^f is very small, the stress-energy and fermionic currents are determined up 
to small correction by the fermionic two point functions. Hence it is sufficient, to determine 
transport coefficients to leading order in 1/iVf, to determine the departure from equilibrium 
of the fermionic two point function. This is precisely the purpose of kinetic theory, so we 
turn to it next. 

It is not our intention in this paper to derive kinetic theory; we feel that the groundwork 
has been adequately established elsewhere [22,23,6]. However, we will outline the parametric 
argument that kinetic theory is applicable to the current problem. The easiest way to see 
that a kinetic theory description is possible is to write down all bubble graphs which appear 
at O(Nf) and at O(iV f ). These are shown in Fig. 1. Cutting one line in a graph gives the 
propagator to the order of interest. The leading order diagram is an undecorated fermion 
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loop. Cutting it gives the bare fermionic propagator, which tells us that the fermionic 
propagator is at leading order a free propagator. (This is obvious since the lowest order 
fermionic self-energy is 0(g 2 ~ 1/iVf).) Cutting the O(iV f ) graph on a gauge boson line 
shows that the gauge boson propagator is at leading order a resummation of one loop self- 
energy insertions, while cutting on a fermionic line shows that the leading fermionic self- 
energy correction involves such a resummed gauge boson line and a bare fermionic line. (The 
gauge boson self-energy is large and requires resummation because Nf different fermions can 
run in the loop, and we must treat g 2 N { as 0(1).) 

This means that it will be possible to use kinetic theory to determine transport coeffi- 
cients in this theory, with errors suppressed by 1/iVf, provided that we carry out a resum- 
mation of the gauge boson self-energy wherever a gauge boson line appears in any process. 
It is incorrect and unnecessary to include gauge bosons as degrees of freedom of the kinetic 
theory; incorrect because they do not have a sharp mass shell, and unnecessary because there 
are so few gauge boson degrees of freedom, so they have a negligible role both in carrying 
conserved quantities and as out states in scattering processes. 

The set of scattering processes needed in the kinetic theory, neglecting 1/iVf corrections, 
is very small, and it will turn out that the scattering processes are simple enough that they 
can be treated without uncontrolled approximations. Note that the large N f expansion is 
much more restrictive than the large N c expansion in QCD; whereas at large N c all planar 
diagrams survive at leading order, at large Nf the structure even of the first subleading 
diagrams is still very simple. In particular it is only at 0(1/ Nf) that diagrams containing 
nonabelian 3-point vertices first appear. This is why, at the order of interest, there is no 
difference between treating QCD and QED. 

Another requirement for a kinetic treatment to be useful is that it is possible to expand 
the space dependence of the two point function in gradients and drop terms with more than 1 
gradient. This is permissible in determining transport coefficients because they by definition 
describe the plasma's response to an arbitrarily slowly varying disturbance. In the case of 
shear viscosity, we consider an arbitrarily slowly varying, divergenceless velocity field. In the 
case of fermionic particle number diffusion, we consider a slowly spatially varying chemical 
potential /i, which we take to be small 3 , C 1. Slowly varying means that the gradients 
must vary on scales well larger than f3Nf. 

In this regime we can write down a Boltzmann equation for the fermionic population 
function (1 particle density matrix) f s (p, x, t) (with s a spin, flavor, particle/anti-particle 
index), describing its time evolution: 



d d d 

777 " — h -t 1 cxt ' 7J 

at ax op 



r(p,x,t) = -0[/], (2.3) 



with v p = p the particle velocity, F ext an external force (if any), and C the collision inte- 
gral, discussed more below. Neither shear viscosity nor fermionic number diffusion require 



3 Nothing prevents considering the case where the chemical potential is large, fi > T, provided <C Nf 
(to avoid complicated physics at the Fermi sphere). However, for \x > T the transport coefficients must be 
reported as a function of g 2 Nf and f3/i, and the energy density becomes /i dependent, so spatial gradients 
in jU also cause bulk flow. We choose {3/j, <C 1 to avoid these complications. 
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including F ext , so we drop it henceforth. The collision integral is zero for the equilibrium 
population function 



/ a (p) = [l + exp(-^K - q s ^)\ , (2.4) 

where q s is the charge of species s under the fermionic number under consideration. The sign 
of q s is opposite between particle and anti-particle. We will expand in the small departure 
from equilibrium, which is justified since the size of the spatial scale of variation of /^(x), 
/z(x) is very large; 4 

f s (p, x, t) = / (p, /3(x), /x(x)) + / (l-/o) 5f s (p, x, t) . (2.5) 

5f s is the quantity we are after, since the transport coefficients are determined by it. The 
nonequilibrium piece of the stress tensor is 

Ta - k 3 V = E/ IPl Vih /o(l-/o) */'(p) , (2-6) 
and the current of species s is 



# = E/ H? s ft/o(l-/o)5f(p), (2.7) 

with the sum over particle and anti-particle. Hence, determining Sf in terms of diUj and 
Vj/i will determine the transport coefficients. 

Since the gradient is taken small, 5f <C fo and we can linearize in 5f. The spatial 
variation, both of fo and 5f, is slow; so we need only keep the gradient terms on the LHS of 
Eq. (2.3) when they act on f . The time derivative term is also irrelevant [20]. The LHS of 
the Boltzmann equation becomes 



0/o(p)[Wo(p)] 



(2.8) 



In writing the term involving ViUj in this form we have used the restriction that we only 
consider divergenceless flow, V • u = 0. The linearized collision operator is rotationally 
invariant, so 5f will have the same angular dependence on p as the left hand side. Since 
the two terms have different angular dependence (the chemical potential involves an £ = 1 
spherical harmonic while the velocity gradient involves an £ = 2 spherical harmonic), the 
two problems can be treated separately. We will define a particle's "charge" q to be q s for 
the case of number diffusion and |p| for the case of shear viscosity. It is also convenient to 
introduce [20] 

r f^-fP'' (diffusion) 

[ V f iPiPj ~ 3°W ' ( shear viscosity) 



4 The distinction between f and 5f is not unique without an additional prescription for Sf and (3^; we use 
the Landau-Lifshitz convention, under which the sum over all particles' Sf gives zero contribution to the 
4-momentum and globally conserved charges. 
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and 



Xt -i( x > — | (ViUj + VjUi - ■ u) . (shear viscosity) ( 2 - 10 ) 



Xi...j parameterizes the source for the departure from equilibrium. The LHS of Eq. (2.8) 
is proportional to X i ...j(x)I i ...j(p). The normalization is chosen so that 7j...j(p)ij....,(p) = 1. 
For two general vectors p and k, 

7 i ... j (p)/ i ... j (k) = P / (p.k), (2.11) 

with Pi the tth Legendre polynomial, and with I the spherical harmonic represented by 
Jj...j, which is the same as the number of indices carries. 

Rotational invariance of the collision integral ensures that we can write the departure 
from equilibrium as 

*/'(p,x) = 3 2 X,... J (x)xl. J (p) , X t,(p) = WpMIpI) • (2.12) 

The function X s (|p|) depends only on the magnitude of the momentum and characterizes 
the departure from equilibrium of the plasma. It has the size of the source for departure 
from equilibrium, JQ...^, and the angular dependence of the departure, scaled out from 
it, and is the most convenient variable in terms of which to write the Boltzmann equation. 
Defining a "normalized source" for the departure from equilibrium, 

Sl. j (p) = -Tq s f (p,x)[l-f (p,x)}I i ... j (f>), (2.13) 

(where again, for diffusion q s is defined above and for shear viscosity q s = |p|), the linearized 
Boltzmann equation, Eq. (2.8), becomes 

5?.. i (p) = (Cxi... i )*(p). (2.14) 

The exact form of the collision operator C will be presented below. The solution Xi-j is 
simply related to the transport coefficients [20]. Introducing the inner product (under which 
the collision operator C is Hermitian and positive semidefinite) 

(f,9)=(3 3 Ej^f a (P)9 a (p), (2.15) 

and considering the case of a fermionic number where all fermions have q 2 = 1, using 
Eq. (2.6) and Eq. (2.7) and the definitions of the transport coefficients leads to 

v = Y^(xij,Cxij) ! 

D=J^r 2 {xi,CXi). (2.16) 

The coefficient in front of D would be 1/3, except that D is defined as —ji/Vijo-, not 
—ji/Vi(J,. The ratio jo/(J> is the charge susceptibility, and for small /i a simple calculation 
gives jo/fi = N { T 2 /3. 
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C. Variational method 



The departure from equilibrium x(|p|) is determined by Eq. (2.14), which is an integral 
equation. (This will become clear when we write the form of the collision integral explicitly, 
Eq. (3.1) below.) Integral equations are difficult to solve unless the integration kernel C has 
a particularly simple form, which is not the case for us. Abandoning an exact solution, we 
can still get a very accurate solution for x, and a more accurate determination of rj and D, 
by expressing the problem as a variational one, writing down a trial function for x(|p|) an d 
solving for a finite number of variational coefficients. The accuracy can then be improved, 
in principle without limit, by enlarging the basis of variational functions considered. This 
technique constitutes a controlled approximation, meaning that it can be made arbitrarily 
accurate. It does not in practice limit the accuracy with which we can extract transport 
coefficients. This strategy has a very long history in the solution of kinetic equations, see for 
instance [24]. Here we follow the approach of [20]. Eq. (2.14) is satisfied at the variational 
extremum of 

Q(X) = (.V,.., , x,..,) - \{xi- 3 , C Xi ...j) , (2.17) 

and the value of the extremum is 

Qmax = 7}(Xi-j » CXi-j) • (2-18) 

This extremal value is simply related to the transport coefficients, presented in Eq. (2.16). 
Since the transport coefficients correspond to the value of Q at a maximum, the error in 
their determination is quadratic in the error between a trial x an d the optimal x- Hence, 
extremizing over a suitably flexible Ansatz for x wm & Ye a ver Y accurate value for the 
transport coefficients. Following [20], we consider the Ansatz 

N 

X s (p)=Y, a mq¥ m) (p), (2.19) 

m=l 

with p — |p| and 0*"^ some set of test functions. Our choice is 

^ (m) ^ = (l^/T^-i » rn = l,...,N, (2.20) 

but other choices also work and lead to the same numerical answer. Note that the (j)^ for 
a small value of N are each linear combinations of the 0*"^ for a larger N, so increasing N 
strictly increases the span of the Ansatz. 

Inserting this Ansatz for x m to Eq. (2.17) turns it into a quadratic equation for the 
coefficients a m ; 

N N 

Q[{ a m}] — ^ Am S m — | ^ dmC mn a n , (2-21) 
m=l m,n=l 

with S m given by S m = (Si...j, and similarly C mn = Cq^ty, where 4™j(p) = 

/ i ... :) (p)0( m )(p). Maximizing Q is now a trivial linear algebra exercise which gives a = C^ 1 S 
and 
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FIG. 2. Scattering processes which must be considered at leading order . 

Qmax = ^ T C^S, (2.22) 

where a = \\a m \\ and S = \\S m \\ are the iV-component coefficient and source vectors, re- 
spectively, in the chosen basis, and C = ||C mn || is the (truncated) collision matrix. The 
determination of transport coefficients then hinges on the accurate evaluations of the in- 
tegrals (Si...j,q<f>\™J) and Cq(f)\ n }j^ . The integral involved in determining S m is 1 
dimensional and is easy to evaluate, quickly and accurately, by numerical quadratures; we 
do not discuss it further. Evaluating C mn is the topic of the next section. However it is 
already possible to determine the power of Nf which will appear in the final answer. The 
quantity S contains a sum over all species but no powers of the coupling constant, so it 
is 0(N { ). The collision integral, since it arises from an 0(1) bubble diagram, must be 
Cmn ~ -^f - The appropriate powers of T follow easily on dimensional grounds. Therefore 
r] oc NfT 3 , and D oc Nf/T (because there is an explicit negative power of Nf in Eq. (2.16)). 

III. COLLISION INTEGRAL 
A. Relevant diagrams 

To make more progress we must analyze what collision processes are important in this 
theory. They are given by cutting the fermionic self-energy diagram of Fig. 1, see Fig. 2. 
Two scattering processes are important, t channel scattering from a particle or anti-particle 
and s channel annihilation and creation of a new pair. The external lines in these processes 
may all be treated as massless, on-shell particles, with uncorrected dispersion relations, since 
the self-energy of the fermion is 0(g 2 ~ 1/iVf). This is not so for the gauge boson propagator 
appearing in the diagram, which must be treated with care. 

Note that only 2 <-> 2 processes appear in the collision integral to the order of interest, 
and that the external lines are all fermions. The set of diagrams required is much simpler 
than the set which are needed at leading order in g when not adopting the large Nf expansion. 
In that case, Compton scattering and annihilation to gauge bosons are already present at 
leading logarithmic order [20], and at leading order interference diagrams between external 
leg assignments of the s and t channel processes are also present. Similarly, in the current 
context, bremsstrahlung emission is also suppressed. These "notable for their absence" 
diagrams are summarized in Fig. 3. In each case there is a clear power counting reason for 
the diagram's absence. Interference between external legs is only possible when all external 
lines are of the same species; but this only occurs in l/N { of scattering events. Similarly, 
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Interference 





and 



Bremsstrahlung 



Annihilation 
(or Compton) 



FIG. 3. Scattering processes which may be important at leading order in the coupling <?, when we do 
not expand in Nf, but which are Nf suppressed and can be neglected here. Dashed lines on the bubble 
graphs show where they are cut to give the diagrams presented. Dotted lines in the scattering processes 
indicate that gauge bosons are not propagating states in the kinetic theory, and we must consider the process 
including the states a gauge boson eventually scatters against or decays into. 



Compton scattering or annihilation to gauge bosons involve a sum over only one, not two, 
fermionic species indices 5 . Another way to see the suppression of these processes is to 
consider what bubble graph must be cut to produce them; this is also shown in Fig. 3, and 
in each case the bubble diagram is 0(Nf 1 ). 



B. Form of C mn 

Since only 2^2 processes are important, the collision integral can be written [20] 



< c *->» = 5 S / wf^^* + k - p '- V)S{ " + * 



P 



k') 



A^( P ,k,p',k') fS(p) jt(k)[l-f^ P ')][l-f^k')] 
X-.,(P)+Xt,(k)-X-.,(P')-Xt,(k') 



(3.1) 



Our notation is that the incoming momenta are p and k and the outgoing momenta are 
p' and k'; and as before, k = |k|, and so forth. A4"^(p, k, p', k') is the matrix element (in 
relativistic normalization) for (ab) to go to (cd). We have already enforced on-shell delta 
functions on the external 4-momenta using free massless dispersion relations, so the energies 



5 In fact, since only the fermions may be treated with kinetic theory, a gauge boson is never allowed as an 
outgoing particle; instead we must draw the eventual scattering the gauge boson undergoes and include those 
fermions as external lines. Since the gauge boson inverse propagator has an 0{g A Nf ~ 1) imaginary part 
where the real part vanishes, the integration over the gauge boson "mass shell" is not dangerous. Drawing the 
diagrams including the fermionic lines makes clearer the absence of Compton and bremsstrahlung processes 
in the current analysis; they really represent either 3 <-> 3 or 2 <-> 4 processes. If we were considering bulk 
viscosity these processes would be important. 



11 



of the external particles equal the magnitudes of their momenta. The leading factor of (1/2) 
corrects for the c <-> d double counting of the sum, or for the outgoing symmetry factor 
when c = d. 

The quantity we want is C mn , that is, Eq. (3.1) with xh replaced with q s (f>\ m }j, multiplied 



by ?°0i"l-(p), integrated over p, and summed over a. The resulting phase space integral and 



matrix element are symmetric on interchange of the (abed) indices, so we can also symmetrize 
the factors, yielding 



abed ' 



x |AOp,k,p',k') / a (p)/o 6 W[l-/o C (pO][l-/o(^)] 

9°*&(p) + ^S(k) - q c ^(p') - (zV&Ml ■ (3-2) 



x 



(p) + ? (k) - <f<^ (p 1 ) - /0ft - 

It is clear from this expression that C is positive semidefinite. For £ ^ 0, it is in fact positive 
definite. 6 

The matrix element gets 3 contributions; the s channel contribution, the t channel con- 
tribution for scattering of like sign (both particle or both anti-particle) fermions, and the t 
channel contribution for scattering of opposite sign fermions (a particle and an anti-particle). 
The matrix element for the s channel process, summed over all species indices (abed) and 
including the symmetry factor of (1/8), is 

s channel : ± ]T \M% * = AtfNtf "( , (3.3) 



abed 



for the U(l) theory. We discuss the nonabelian generalization in the next section. The 
meaning of "[t 2 + u 2 /s 2 }" is a quantity which equals this (with s,t,u the usual Mandelstam 
variables) in the small g 2 N{ limit, and is given in complete detail in Appendix B Eq. (B7). 
The expression for t channel, like sign scattering is 



i " / 2 i 2 \ " 



^ abed 



t channel like sign : - £ \M% = A(g 2 N { ) 2 [ ) . (3.4) 



The expression for opposite sign scattering is identical. Again, the quantity "[s 2 + u 2 /t 2 ]" 
equals this in the small g 2 Nf limit, but is given by a more complicated expression presented 
in Appendix B, Eq. (B8). 



C. Doing the integrals 

To compute transport coefficients it remains to perform the integrals presented in 
Eq. (3.2). There is a 12 dimensional integral to perform. Four of the integrals will be 



6 It is not positive definite for I = 0, in which case the collision integral vanishes if <j>(p) is a constant, until 
we include higher order, number changing processes. This is discussed more in [6], and is the reason that 
evaluating bulk viscosity is difficult. 
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performed by the energy-momentum conserving delta functions, and three are trivial be- 
cause the way we structured Eq. (2.8) makes Eq. (3.2) invariant under global rotations. 
This leaves 5 integrals. Appendix A shows two nice ways of parametrizing these remain- 
ing integrals. One way, due to [13], is convenient for t channel exchange, and another is 
convenient for s channel exchange. 

We can write Eq. (3.2) as a sum of two 5 dimensional integrals, one containing the s 
channel contribution to \ Ai\ 2 and the other containing the t channel contributions, by using 
Eq. (A5) and Eq. (A14). The matrix elements are given by Eq. (3.3) with Eq. (B7), and 
twice Eq. (3.4) with Eq. (B8). The product of functions is evaluated for shear viscosity 
using 

|p|0!7 ) (p)|k|0j ) (k)=pA;0( m )(p)0W(A;)P 2 (cosV), (3-5) 

and similarly for the other terms when the product of 0's is expanded. Here Pii^) — 
(3x 2 — l)/2 is the second Legendre polynomial, and all the cosines of angles are given 
in Appendix A. For both the s and t channel part of the scattering integral, there is an 
integration over an azimuthal angle which can be done analytically, leaving a 4 dimensional 
integral, over the frequency uj and momentum q carried by the gauge boson line and two 
external momenta. The integrand is long and unenlightening and can be obtained from 
the equations listed above. We perform it by numerical quadratures using adaptive mesh 
refinement methods. The accuracy with which the integrals can be performed limits the 
accuracy of the determination of the transport coefficients, but it is not too difficult to 
achieve 0.1% accuracy in the finally determined transport coefficients, for the full range of 
the coupling constant g 2 N^ that we consider. Note that, besides the explicit g 2 N^ dependence 
of Eq. (3.3) and Eq. (3.4), there is also quite complicated g 2 N{ dependence in the matrix 
elements because of the gauge boson self-energies. Therefore the full integration must be 
re-performed for every value of g 2 Nf under study. Note also that g 2 Nf requires (vacuum) 
renormalization. We regularize using the MS scheme. We also choose the renormalization 
point 7 p, to be the "dimensional reduction" value [25] /xdr = 7rTexp(— 7 E ), with 7 E = .577. . . 
the Euler-Mascheroni constant. This choice is arbitrary and represents a choice in the 
presentation of the results only. Its motivation is that, at this renormalization point only, 
the free energy of an infrared electromagnetic field is given by / d 3 x(E 2 + B 2 )/2. 

We also cannot treat arbitrarily large g 2 Ni because the presence of the Landau pole 
in the theory becomes problematic. Requiring that the Landau pole occur further into 
the ultraviolet than [g^g^l = (40T) 2 requires g^/W{ < 4.4. A Landau pole this far in 
the ultraviolet does not affect our calculation because the contribution to Eq. (3.2) from 
integration regions with \u 2 — q 2 \ large enough to encounter the Landau pole is exponentially 

small. In fact such regions were never sampled by the adaptive integration routine. For 

l li 

gN { = 3, the energy scale where the pole appears is already above 1000T. 

For the case of fermion number diffusion, there are some additional cancellations which 
make the structure of the integrand slightly simpler than for shear viscosity. Unlike shear 
viscosity, the charge q s for number diffusion is opposite for a particle and its anti-particle. 
Further, the t channel matrix element to scatter from a particle equals the matrix element 



7 Do not confuse the MS renormalization point \i with the chemical potential \i. 
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to scatter from its anti-particle (up to 1/iVf suppressed corrections). Hence, for t channel 
exchange, summing on like and opposite sign scatterings cancels all "cross-terms" in Eq. (3.2) 
with bothp and k line arguments for 0's. Using also the p <-> k symmetry of the integration, 
and the fact that the incoming and outgoing particles a and c are the same type and have 
the same charge, one may substitute 



\M\ 2 {t channel) 

abed 



, (m) 



P ) + gVS m) (k)-g 

! n) (p) + gV! n) 

(m) / \ ±i m ) 



X 



^2|7W| 2 (t channel)(g a ) 

abed 



?>r(p)+?vr(k)-? c 



(p) 



(p') 



t\p' 


)- q d < 


Pi 


i n \p) 


-q d 4 


!"*( 


S n) ( P ) 


- <P \ 


(pO 



(k' 
k')" 



For the s channel exchange a similar simplification occurs because in that case q c 
the integral is symmetric on p' <-> k' (as well as on (p, k) <-> (p', k')); so 

£ \M\ 2 (s channel) W m) (p) + ^(k) - <f0fV) - rt M (k') 



(3.6) 



-q d but 



abed 



q^ip) + ?VS n) (k) - ^(p') - g d 0l n) (k') 



^2\M\ 2 (s channel) (q 

abed 



a\2 



p)-0, (m) (k) 



p)-0S n) (k) 



(3.7) 



(The same symmetries were used in [20] to allow a simple leading log calculation of the diffu- 
sion constant.) For both the s and t channel integrations, there is an azimuthal integration 
which is trival, and the two external momentum integrations factor and can be performed 
separately (rather than being nested). This simplifies somewhat the numerical integrations. 

Further, the cancellations just described show that, for number diffusion, the way that 
the departure from equilibrium in one species evolves decouples from the departure from 
equilibrium in every other species. 8 This is why all the fermionic diffusion constants are the 
same, unless the total fermionic number is coupled to a gauged U(l) field. 

IV. ALL-ORDERS RESULTS 



The last two sections describe how the transport coefficients, at leading order in 1/iVf but 
all orders in g 2 Nf, may be computed for U(l) gauge theory. The results appear in Figures 
4 and 5. Here we discuss the results and generalize them to QCD, with N c Nf so as not 
to spoil the expansion we have used. 

First we remark that, while there is no "diffusion constant" for the total fermionic number 
density when it is coupled to a U(l) gauge field, the would-be diffusion constant instead 
determines the electric conductivity for the U(l) gauge field [20]: 



(4.1) 



s Note that this would not happen if we considered chemical potentials [i > T. 
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FIG. 4. Shear viscosity, computed to all orders in g 2 Nf. The renormalization point /i in the MS scheme 
is p, = 7re~ 7E T, the "dimensional reduction" answer. 
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FIG. 5. Diffusion constant, computed to all orders in g 2 Nf . The renormalization point is the same as in 
Fig. 4. 
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Therefore the results in Fig. 5 also constitute a result for the electrical conductivity of this 
theory. For this interpretation the choice of renormalization constant is essential, since only 
for p, = /i DR does the free energy density stored in an electric field equal E 2 /2. 

In presenting the transport coefficients in Figs 4 and 5, we factor out the leading g~ A 
behavior from the transport coefficients; but there remains nontrivial g 2 N dependence. This 
is because the size of the gauge boson self-energy, relative to the tree gauge propagator, 
depends on g 2 N{. It is not surprising that this is important for g 2 N > 1. In fact the 
dependence continues down to arbitrarily small g 2 Nf, where the transport coefficient scales 
as 1/ ln(l/g 2 iVf). This is because the t channel contribution to the scattering integral, 
Eq. (3.2), has a logarithmic small q divergence which is cut off by the presence of a self- 
energy on the gauge boson line. This is discussed further in the next section. 

The results presented in Fig. 4 and Fig. 5 are for U(l) gauge theory. Consider instead 
a general nonabelian theory, with Nf flavors of fermions in a common representation of 
the group, of dimension d F and with quadratic Casimir C F , and trace normalization T F = 
C F d F /d A , with d A the dimension of the adjoint representation. For SU(iV c ) gauge theory 
with fundamental representation fermions, C F = (N 2 — 1)/2N C , d F = N c , and T F = 1/2. 
The total number of fermionic degrees of freedom is then 4d F N{ , rather than 4Nf in the U(l) 
theory. (The 4 is because we consider Dirac fermions; there are two spins, and particle/anti- 
particle.) Therefore, the sum in Eq. (2.15) includes 4Nfd F terms, rather than 4iVf terms as 
in the U(l) theory, so S m is multiplied by d F . Meanwhile, the gauge field self-energy has 
the replacement 

2 9 2 N f T 2 2 g 2 (T F N i + T A )T 2 g 2 T F N { T 2 
™d = — 3 — m D = 3 3 • ( 4 -2) 

The expression with Ta (Ta = Ca the adjoint Casimir, Ta = N c in SU(A r c ) gauge theory) 
holds generally to leading order in g 2 in QCD, with the Ta term arising from gauge boson 
loops; in the last expression we have enforced the large Nf limit to eliminate this nonabelian 
term. The group theory coefficient for C mn is (d F C F /g?a) • Also, the fermion number suscep- 
tibility, which goes into the determination of D, is multiplied by d F ; hence D is smaller by 
this factor. Therefore, to use the figures to get the transport coefficients for a nonabelian 
theory, one should interpret the x axis in Fig. 4 as 

gN 1 , 12 => g{T F Nffl 2 , (4.3) 

and replace the y axis label with 

9 4 9 A Cl „ 4Ar ^ „ g 4 NfC F d F T 

In SU(A^ C ) QCD with fundamental fermions, T F = 1/2, d F = N c , and C F jd^ = {N 2 — 
1)/(4:N 2 ). In converting D into a using Eq. (4.1), the g 2 appearing in that equation should 
be the U(l) gauge coupling, and there is an additional factor of d F , 

e 2 d F N(T 2 

a = D , (4.5) 

with e the fermionic coupling to the U(l) gauge field. 
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V. LEADING ORDER CALCULATION 



One use of our results is to allow a comparison against a leading order calculation which 
relies on the hard thermal loop expansion. To date, only leading logarithmic results are 
available. Repeating the analysis of [20] in this context gives at leading logarithmic order 

91.1725 1 _ 61.0727T 3 1 

~ Nfg 4 T X ln(T/m D ) + 0(1) ' V ~ g 1 X ln(T/m D ) + 0(1) ' ( ' 

To study better the reliability of the HTL expansion, we present here a leading order calcu- 
lation. We perform it in the spirit of other leading order HTL calculations, which is to say, 
all self-energies are replaced with just the thermal part of the self-energy in the HTL approx- 
imation, wherever that treatment is parametrically justified. Note that the HTL vertices 
all vanish at leading order in 1/Nf and are not needed here. In this section and this sec- 
tion alone, we will treat g 2 N{ as a small quantity which can be expanded in parametrically. 
Equivalently we treat the Debye mass, given by = g 2 NfT 2 /3, as m D <C T. 

It turns out that both the t and s channel exchange parts of the collision integral depend 
in a nontrivial way on the presence of a self-energy on the gauge boson line. We will discuss 
each in turn. 



A. Leading order t channel exchange 

It has long been appreciated that there is a subtlety in the t channel exchange process, 
and that inclusion of a self-energy on the gauge boson lines is necessary, not only to get the 
correct leading order behavior, but to get a behavior which does not suffer from infrared 
divergences in the evaluation of the collision integral [13]. This is because of the small 
denominator provided by the propagator when p — p' is small. As shown in Appendix A, 
the t channel phase space integrals can be written [13] (Eq. (A5)) 

0& /.qo ,.q /.qo />oo i>2n 

' dq du dp dk d<p , (5.2) 



2 9 7T 6 Jo J-q Ji=f J^L JO 

with p and k the energies of the two incoming particles, and oo and q the energy and 
momentum carried by the exchanged gauge boson. When the particle energies are O(T) but 
the exchange momentum is mD < g C T, we can expand in small q/p and simultaneously 
drop the self-energy effects on the propagators, so that Eq. (B8) becomes 

a f s 2 + u 2 \" s 2 + u 2 8p 2 k 2 (l-cos6) 2 

At first sight the integral over q is severely infrared divergent. However, two more powers of 
q arise from cancellations in the expressions of the form 0j...j(p) — 0j...j(p') ... in Eq. (3.2) 
[20], because p— p' = q is 0(g) and <pi...j is a smooth function of its argument; so the 4>i-jip) 
and 0j...j(p') terms cancel up to an 0(q) remainder. Therefore, the q integration is really of 
the form 

dq / duo— , (5.4) 
J—q q 
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which is logarithmically sensitive to the gT scale. The upper end of the logarithm is cut off 
because the expansion q <C p is no good there. The lower end of the integration range is 
only cut off by the inclusion of the self-energies in the scattering processes. Therefore it is 
necessary to include those self-energies to get a finite answer. However, for m D <C T, they 
only play a role when q < tod <C T, which is clear from the form of the matrix element, 
given in Eq. (B8). 

At leading logarithmic order it is only necessary to determine the coefficient of the log 
divergence arising from the C g < T momentum range. Such a calculation was carried 
out in [20], where it was also not necessary to appeal to a large N { expansion. However, a 
leading order calculation requires a treatment of both the q ~ and the q ~ T momentum 
ranges, as well as the intermediate range. For such a leading order calculation we must 
perform the full integration presented in Eq. (3.2). However we may simplify the treatment 
by replacing the full self-energy on the gauge boson line with its HTL counterpart. That 
is, we may drop the vacuum contribution to the self-energies, and make the simplifying 
approximations for the thermal parts which are presented in Eq. (B12). We may do this 
throughout the integration range of Eq. (3.2). The vacuum part is uniformly suppressed by 
g 2 N{ (up to logarithms). At small transfer momenta q ~ mo, the corrections to the HTL 
approximation for the thermal contribution to self-energies is Oivn^jT 2 ~ g 2 Nf), and the 
approximation is valid. At large momenta q ~ T, the HTL approximation is invalid, but 
the whole thermal self-energy is an 0(g 2 N { ) correction and is negligible anyway. 

Our method for achieving a leading order determination of the integrals in Eq. (3.2) 
is therefore to replace the self-energy on the gauge boson line with the HTL self-energy, 
over the full momentum integration range. The integral is then performed without further 
approximation (such as expansion in q <C p). 

We have also tested the sensitivity to a slight modification of the matrix element, moving 
the cos(/> independent piece in Eq. (B8) which is multiplied by (q 2 -l-Hr) -2 over to the piece 
multiplied by (q 2 + Hl) -2 , so appearances of the transverse self-energy coincide with powers 
of cos(/>. This modification makes an 0(g 2 N{) change because the piece moved is 0(g 2 N{) 
suppressed where the difference in self-energies matters. It coincides with the separation 
between transverse and longitudinal self-energies presented in [14], and is the easiest thing 
to implement, especially if scatterings with gauge bosons are important (in a nonabelian 
theory and not at large iVf). We find the difference is relatively small, comparable to the 
difference between the leading order at the dimensional reduction renormalization point and 
the all orders answer. This is not the procedure we used to get the results presented in the 
figures to follow. 

B. Leading order s channel exchange 

Naively, the s channel exchange contribution to Eq. (3.2) should not need any self- 
energy insertion at all on the gauge boson propagator. After all, the exchange energy is only 
small when both incoming fermions have small momentum, which is a strongly phase space 
suppressed region. Since we are working with fermions there is also no enhancement from 
the statistical factors in this region. However, while the self-energy is in fact negligible at 
generic momenta, this is not so when the incoming fermions are hard but collinear, p,k ~ T 
but 1— p • k ~ g 2 N { , or equivalently u> — q ~ g 2 NfT. Furthermore, this region turns out to 



18 



FIG. 6. Diagram for production of an on-shell gauge boson. 



contribute 0(1) of the contribution of the entire s channel piece of Eq. (3.2). The reason 
is that the thermal self-energy lifts the transverse gauge boson mass shell from lying on the 
light cone u> 2 — q 2 = to lying in the timelike region, lo 2 — q 2 ~ m^. Since the fermions still 
propagate at the speed of light, it is kinematically allowed for a fermion pair to annihilate 
to an on-shell gauge boson, which propagates for some time and then breaks up back into 
fermions; that is, there is resonant enhancement of the scattering cross section on the gauge 
boson "mass shell". The phase space suppression is compensated for by the on-shell near 
singularity from the propagator. 

At finite g 2 N{, this does not lead to any divergence in the the matrix element at any 
value of lo ^ q, because the gauge boson self-energy has a nonzero imaginary part everywhere 
except the light cone. The imaginary part "smears out" the gauge boson mass shell from a 
delta function into a sharply peaked Lorentzian. We now estimate the matrix element near 
resonance. The relevant part of the matrix element, Eq. (3.3), is the term involving the 
square of the transverse propagator, which is of the form (see Eq. (B7)) 



/ 2 2\2 

\M\ 2 (s channel) ~ f- ~ 9 >- x 0(1) . (5.5) 
\lu z — q — 11t| 

The real part of IT T is ~ m^, which is relevant for (lo — q) ~ m^/(a; + q), that is, within 
g 2 T of the light cone. At timelike momentum but near the light cone, the leading order 
expansion for the real and imaginary parts of IIt, keeping the largest real and the largest 
imaginary part, is 

Ht(u, q) < — ~ — ^ Jj2 ~ ^ ^ ~~ ^ nerma ^ P ar t] • (5.6) 

The thermal imaginary part is of the opposite sign as the vacuum imaginary part, and is 
strictly smaller. Its form is rather complicated even in the small (ou—q) limit. Note that the 
imaginary part vanishes as the momentum approaches the light cone, q — > lo. The imaginary 
part vanishes in the standard HTL approximation to the self-energy. 
The location of the mass shell, for lo ~ T ^> m D , is 

(5.7) 

just as if the gauge boson had a mass squared of The half width of the mass shell is 

g 2 N { ml 

Aqi/2 ~ ' (5 ' 8) 
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which is very narrow, 0(g A N 2 T). This gives a phase space suppression to the resonance 
region. The value of the matrix element on resonance is 

which is exactly enough to make up for the phase space suppression and make the contri- 
bution to the scattering rate, from exchange of an on-shell gauge boson, comparable to the 
integral over all exchanges of off-shell gauge bosons. 

An easier way of seeing this result is to consider the process of on-shell gauge boson 
production, shown in Fig. 6. This process contains an explicit g 2 Nf, and it is also phase 
space suppressed for kinematic reasons, by O^m^/T 2 ), which is 0(g 2 Nf). Therefore it 
contributes 0(g i Nf) to Eq. (3.2), which is parametrically of the same order as the other 
scattering processes, even when the coupling is taken small so the gauge boson propagates 
close to the light cone. 

We parenthetically remark that the "formation timescale" for the collinear annihilation 
process just described is ~ g/m^ ~ 1/ ' g 2 N{T. This is shorter than the mean time between 
scatterings for the fermions, ~ 1/ 'g 2 T, by a factor of 1/N { <C 1. Therefore the annihilation 
process, although it involves collinear physics, is not sensitive to scatterings of the fermions 
at leading order; there is no Landau-Pomeranchuk-Migdal suppression. However, if we were 
speaking of QCD or QED outside the large Nf expansion, this would not be the case. The 
collinear pair annihilation to a gauge boson process would occur at a rate of 0(1) importance 
to transport coefficients, but the rate itself would get 0(1) corrections from scatterings on 
the incoming legs. This has recently been discussed, in the context of photon emission from 
the quark-gluon plasma, by Gelis et. al. [26]. 

There are two strategies to compute the leading order in g 2 N{ contribution of s channel 
scattering to Eq. (3.2). One method is to compute the s channel scattering diagram, with no 
self-energy insertion on the gauge boson propagator, and to add to it the rate of the // <-> 7 
process just discussed, computed using the corrected gauge boson dispersion relations and 
treating the gauge boson as an external line. This method means that we are allowing gauge 
bosons to be external states in scattering processes, which means that they are being included 
as kinetic degrees of freedom. Such a kinetic description of the fermions and gauge bosons 
is possible at leading order in g 2 N { and should be distinguished from the kinetic treatment 
of the fermions alone, which is accurate to 0(1 /Nf) corrections as discussed above. 

A kinetic treatment of both fermions and gauge bosons requires considering the departure 
from equilibrium of the gauge bosons. For number diffusion this is trivial; the departure 
vanishes, as the gauge bosons carry no fermionic number. For shear viscosity it is not trivial, 
though. Therefore we turn to an alternative, which gives the same result. It is to compute 
the rate of the // <-> // process at finite g 2 N { , using either the full self-energy or the 
approximation of Eq. (5.6) above, and to extract numerically the small g 2 Nf limit. The 
convergence to a small g 2 N{ limit is very good, with power suppressed corrections. 

We emphasize, however, that it is not correct to compute the leading order scattering 
contributions from s channel exchange by simply dropping the self-energy on the gauge bo- 
son line. This misses the contribution from resonant scattering through an on-shell gauge 
boson. It is also not correct to insert the HTL self-energy on the gauge boson line. The HTL 
approximation to the self-energies is correct at leading order, not only at soft momentum 
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q, uj ~ rriB, but at nearly lightlike momentum 9 u 2 — q 2 ~ m 2 ^, and this is the important 
regime. However, in the HTL approximation the self-energy has no imaginary part at time- 
like momentum. The presence of a small imaginary part is essential to prevent a divergence 
in the scattering rate, so it is necessary here to go beyond the HTL approximation and get 
the leading nonzero imaginary part. 

Numerically, the contribution of the on-resonance region is comparable to the off- 
resonance s channel region for the case of number diffusion, and is a few times smaller 
for the case of shear viscosity. The importance of all s channel scattering processes, relative 
to t channel contributions, is numerically small for both transport coefficients. Throughout 
the range of coupling g we have considered, the t channel scattering part of Eq. (3.2) domi- 
nates the s channel part by at least a factor of 3 for number diffusion and a factor of 20 for 
shear viscosity. The domination grows larger for smaller values of g. The dominant role of t 
channel scattering processes in setting the shear viscosity is also found at leading log order 
[20]. 

C. Leading order results 

We compare the leading order calculation to the complete calculation in Fig. 7 and Fig. 8. 
The errors due to imprecision of numerical integrations are smaller than the line thicknesses 
in the figures. Note however that there is a formally 0(g 2 Nf) ambiguity in the leading order 
calculation. The self-energy is taken to be only the HTL thermal one, without the vacuum 
contribution. Therefore, in the leading order calculation, the value of g 2 does not vary with 
q^q^ = q 2 — u 2 , but is fixed. In the complete calculation the coupling is the running coupling. 
To compare the results of the leading order and complete calculations we must decide at 
what renormalization point to evaluate the running coupling constant, to insert as the fixed 
g 2 of the leading order calculation. This is a normal problem in gauge theories; a calculation 
made to some order in a coupling constant carries renormalization point dependence at one 
higher order in the coupling. 

There is a particularly natural choice for the renormalization point, which is the one for 
which tree expressions give the correct infrared, thermodynamic description of the gauge 
field physics. This is the "dimensional reduction" value /2 = ne~" m T [25]. We have chosen 
this as a "central" value for the presentation of results, and we illustrate the dependence on 
Jx by varying renormalization point by a factor of e on either side. 

An alternative approach is to settle the renormalization point dependence by including 
the real part of the (very simple) vacuum self-energy contribution everywhere, in addition 
to the HTL approximation of the thermal contribution. This is equivalent to using g 2 (j± = 
\u 2 — q 2 \) as the value of g 2 when the exchange momentum is (u;,g)-that is, allowing the 
gauge coupling to run as the exchange momentum varies inside the integral in Eq. (3.2). We 
therefore term this a "renormalization group improved" leading order calculation. It is not 
correct beyond leading order in g 2 Nf, but it is free of renormalization point ambiguities. It is 
also shown in Figs. 7 and 8. The RG improved result systematically lies above the p, = /2dr 



9 We thank Tony Rebhan for pointing this out to us. 
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Shear viscosity, N f >>l 




Leading Log (2 consts) - 

Q I I I I I I I I I I I I I I I I I I I I I I 

12 3 4 

N? /2 g(M DR ) = (m D /T)V3 

FIG. 7. Shear viscosity computed with various approximations. The leading order calculation has an 
0(g 2 N{) ambiguity associated with the renormalization point chosen for g 2 ; we show the result for three 
choices of renormalization point, p, = /Idr, and e and 1/e times this, and "RG improved." The leading log 
result is shown for two constants under the log. The curve which grossly disagrees with all others is using 
0(1) = in Eq. (5.1), the other leading-log curve picks the 0(1) constant to reproduce the leading order 
result at the smallest value of gy^W- 
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Diffusion Constant, N f >>l 




Leading Log (2 consts) - 

i i i i I i i i i I i i i i I i i i i I i 

12 3 4 

N, 1/2 g(M DR ) = (m D /T)V3 

FIG. 8. Same as Fig. 7 but for flavor diffusion. 
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value, and turns out to be a worse approximation over the entire range considered. Note 
the strikingly good performance of the p = /xdr renormalization point, leading-order result, 
especially at small coupling. We interpret this as the result of "getting the renormalization 
point right" ; the leading thermal correction, left out when we make the HTL approximation, 
cancels the vacuum running of the coupling and effectively sets g 2 ~ (? 2 (/Udr) in the soft 
exchange momentum region. 

The figures also show the result of the leading-log calculation. There is an even more 
severe ambiguity in the leading-log calculation than in the leading order one; we do not 
know the 0(1) constant in Eq. (5.1). We display two values. The curve in each figure 
which lies above all other curves is setting 0(1) = in Eq. (5.1), while the curve which 
tracks very closely the Ji = /iDRe -1 leading order curve is the result if we choose the 0(1) 
constant to match the leading order calculation at the smallest value of itld/T shown. The 
ambiguity in the leading-log result is large even where the coupling is weak, and the naive 
choice 0(1) = in Eq. (5.1) gives horrible results. However, what amounts to a "next to 
leading log" calculation works surprisingly well, comparably to leading order. 

VI. CONCLUSIONS 

Computing transport coefficients, or any other quantity involving infrared limits of corre- 
lation functions, is a difficult problem in the context of relativistic gauge theories. We have 
shown that, in a rather special limit, this problem simplifies and the calculation, though 
still difficult, is tractable. This limit is the limit of a large number of fermions, Nf >> 1, 
expanding to leading nontrivial order in 1/iVf. In this limit, nonabelian effects become unim- 
portant. Furthermore, the fermions enjoy nearly free dispersion relations, with corrections 
(both dispersion and scattering corrections) suppressed by 1/iVf. This allows a kinetic theory 
treatment of transport coefficients, which we have presented. 

This kinetic treatment is actually simpler than a kinetic treatment of QED or QCD 
would be, if we intended to compute to leading order in g. In particular, in large Nf 
QED, there are no contributions at leading order from Compton scattering and annihilation 
to gauge bosons, or from interference diagrams or bremsstrahlung processes. Nonabelian 
effects are also irrelevant at leading order, so the QED and QCD calculations differ only by 
group factors. The structure of the scattering integral is simple enough that it is possible 
to carry the calculation out nonperturbatively in g 2 N { , though this demands some care in 
the treatment of the scattering matrix element. In particular, the hard thermal loop (HTL) 
approximation for the thermal self-energy is inadequate, and a general expression for the 
thermal self-energy is needed. 

It is also possible to perform a leading order in g 2 Nf calculation by making the hard 
thermal loop (HTL) approximation. This permits an interesting test of the quality of the 
HTL approximation. Because the physics of this model is missing some complications present 
in QED and QCD (such as bremsstrahlung), we expect the HTL approximation to work, 
if anything, better in this model than in realistic QED or QCD. In fact, corrections to the 
HTL leading order calculation which we made here are parametrically 0(g 2 Nf), while in full 
QED or QCD we believe that gauge boson emission processes and nonabelian effects (QCD 
only) will introduce 0(g) corrections. 
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Our results for the comparison between the complete solution and the HTL approxima- 
tion are presented in Fig. 7 and Fig. 8, which constitute the main results of this paper. They 
show that, for g 2 N { < 1, the leading order treatment is accurate to a few percent. This is 
good news for QED, where the coupling is indeed small. For instance, between the electron 
and muon mass scales, the QED coupling is g^/N{ ~ .30. At this coupling the difference 
between the HTL and full treatments is less than 0.3%. Even at this coupling, however, a 
calculation to leading order in the logarithm of the coupling is of very little value. 

On the other hand, as the coupling grows larger, the leading order treatment becomes 
much less reliable and more renormalization point dependent. In the theory we have looked 
at, the size of the renormalization point dependence of the leading order calculation and 
the size of its error (its difference from the exact result) are comparable. However, the 
"RG improved" calculation shows that the error is probably not due to a poor choice of 
renormalization point, but is a failure of the HTL approximation to correctly represent the 
thermal part of the self-energy. This is mixed news for QCD. 

We can estimate how accurate a leading order calculation of transport coefficients in 
QCD would be by seeing how accurate the large Nf results are at a value of the Debye mass 
mji/T equal to the relevant QCD value. First consider QCD for T at the electroweak scale, 
T ~ lOOGeV. Taking a s — 0.11, iVf = 6 and using Eq. (4.2), we find this corresponds to 
g^/Nf — 2.9. For this value, the renormalization point uncertainty in the leading order value 
is about 30%. This is not so bad, though in full QCD the actual error will likely be larger 
than this. 

The correction is much more severe at T ~lGeV, where the coupling is much larger. 
This value of T gives a g roughly equal to the right hand edge of the figures, where the 
renormalization point sensitivity of the leading order calculation is at least a factor of 2. 
Therefore, at temperatures which may realistically be obtained in heavy ion collisions in the 
foreseeable future, a full leading order calculation is at best a factor of 2 estimate. 
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APPENDIX A: INTEGRATION VARIABLES 

Both s and t channel scattering processes must be computed in the main text. Different 
simplifications of the phase space integration in Eq. (3.2) prove most convenient for s and t 
channel scattering processes. It is convenient to have, as integration variables, the magni- 
tudes of the spatial and temporal components of the momentum carried by the gauge boson, 
since it is by far easiest to express \Ai\ 2 in terms of these. Here we present both choices for 
phase space parameterization we have found useful. Note that we always take the external 
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momenta to be strictly lightlike (massless dispersion relations), meaning po — |p| = P, which 
is justified at leading order in l/N f . 

1. t channel exchange 

First consider the case were p' — p is the exchange momentum. Then in the collision 
integral , Eq. (3.2), it is convenient to use the spatial S function to perform the k' integration, 
and to shift the p' integration into an integration over p'— p = q. We may write the angular 
integrals in spherical coordinates with q as the z axis and choose the x axis so p lies in the 
x-z plane. This yields 
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roo rl rl'K | 

I q 2 dqp 2 dpk 2 dk d cos 9 pq d cos 9 kq / defy— —5(p+k-p'-k') , (Al) 

o J-i Jo pkp' k' 



where q — |q|, and p' and k! become dependent variables, p' = |q + p| and k' = |k — q|. is 
the azimuthal angle of k (and k') [i.e., the angle between the p-q plane and the k-q plane], 
and 9 pq is the plasma frame angle between p and q, cos 9 pq = p • q, etc. 

Following Baym et al. [13], it is convenient to introduce a dummy integration variable 
cu, defined by a 5 function to equal the energy transfer p' — p, so that 

/'OO 
duo 5{uo+p- p') 5{uo-k + k') . (A2) 
-oo 

Evaluating p' — |p + q| in terms of p, q, and cos9 pq , and defining t = uj 2 — q 2 (which is the 
usual Mandelstam variable), one finds 

5(uj + p-p') = ^-5 ( cos 9 pq - - - 0(u; + p) , (A3) 
pq \ q 2pqJ 

S(u-k + k!) = 5(cos6 kq - ^ + " ^) • (A4) 

Here is the step function. The cos 9 integrals may now be trivially performed and yield 1 
provided p > \{q — u), k > \{q + u>), and < q; otherwise the argument of a 5 function 
has no zero for any | cos^l < 1. The integration range becomes 



05 rOC rq 



rOO rq rOO rOQ r2ir 

2%«i'"/./ a i i? *4' it /„*- < A5) 

with p' = p+uj, k! = k—uj. For evaluating the final factor of Eq. (3.2), note that 

I i ... j (p)I i ... j (k) = P e (cos9 pk ), (A6) 

where is the £'th Legendre polynomial, and I = 1 if there is one index (diffusion) and 
£ = 2 if there are two indices (shear viscosity). We will therefore need expressions for the 
angles between all species, as well as the remaining Mandelstam variables s and u, which 
may appear in j\4 2 . They are 
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and 



-t 
S ~2f 
u = —t — s , 



(p + p')(k + k') + q 2 - cos yj (App' + t) {Akk! + t) 



COS 9 pq = 



ijj t 

q 2pq ' 



cos 9 kq 






Q 


2kq ' 


cos 6 pp i 


= 1 + 


t 




2pp' ' 


cos 9 pk > 


= 1 + 


u 


2pk' ' 


cos 6 pk 




s 


= 1 - 


2pk ' 





COS 6p/ q 



COS ^fc/, 



q 2p'q ' 

a; t 
- + 



cos 9 kk > = 1 + 



2kk' ' 
u 



cosfl^ =1 + — , 
cos 6 p >k> = 1 



2^' 



(A7) 
(A8) 



(A9) 



The integration can be done analytically in all cases in this paper, but the other integrals 
cannot. 



2. s channel exchange 

When the exchange gauge boson carries momentum p + k, it is convenient to use a 
different parameterization of phase space. Again the spatial 8 function in Eq. (3.2) is used to 
perform the k' integration, but now the k integration is shifted to an integral over q = p + k, 
the total incoming spatial momentum. Again we use spherical coordinates, with q on the z 
axis and p in the x-z plane. In these variables, the integration measure becomes 

? I SX' (2r)V(p + k " p ' " k ' )Sip + k - p '~ k>) 



poo /"i I 

„ 6 J o <l 2d( l P 2d P V dp' J d cos 9 pq d cos 9 p , q d(j> - 5 (p+k-p'-k 1 ) , (A10) 

where k = |q— p| and k! = |q— p'|, and <fi is the azimuthal angle of p', meaning the angle 
between the q, p plane and the q, p' plane. 

Next, we re- write the energy delta function in terms of u, the total energy; 

PCX) 

S(p + k-p' - k') = du5{u-p-k) 5{u - p' - k') . (All) 
Jo 

In terms of the remaining integration variables, and defining s = uo 2 — q 2 (which again is the 
Mandelstam variable), one finds 

5(uj-p-k) = — 5[cos9 pq - - + ^-)e(w-p), (A12) 
pq \ q 2pq J 

S (u - p '-k') = ^-5 (cos9 p/q -- + -^-) Q(lu - p') , (A13) 
pq \ q 2p'q J 
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which performs the integrals over cos6 pq and cos 6 p t q . There is only a solution provided 
q < u), \2p — uj\ < q, and \2p' — uj\ < q. The integration measure becomes 



P 3 



dw / dq 



u — q 
2 



dp 



LJ — q 
2 



2tt 



2% 6 Jo 

with k = Ld—p, k' = uj—p'. The Mandelstam variables are 
(p - k)(p' - k!) - q 2 



dp / d<p , 



u 



2q 2 
-s-t. 



+ cos (j)\J (Apk — s)(Ap'k' — s) 



The angles with respect to q are 

uj s 



cos 9 pq 
cos 8 kq 



q 

UJ 



2pq 

s 

2kq 



COS 9pi q 
COS 6 k 'q 



LU 

1 

UJ 



2p'q ' 

s 



2k'q ' 



(A14) 

(A15) 
(A16) 



(A17) 



and the remaining angles are as in Eq. (A9). 

Again the integral can be performed analytically, but the other integrals cannot. 



APPENDIX B: MATRIX ELEMENTS 

In this appendix we present the matrix elements for s and t channel gauge boson ex- 
change. These must be calculated without making any expansion in g C p or in small 
self-energy; for our purpose the self-energy must be taken to be an 0(1) correction and all 
momenta must be considered to be the same order. 

The quantity u [t 2 + w 2 /s 2 ]", introduced as a shorthand in Eq. (3.3), technically means 

1 -D^DtfTr(p<r M)^Ai>'l a ■ (Bl) 

Following the previous appendix we denote the plasma frame frequency and momentum 
carried by the gauge boson propagator as uj and q. The retarded gauge boson propagator 
D T ^ is most conveniently expressed in Coulomb gauge-the result for Eq. (Bl) being gauge 
invariant: 

^•4 = 7=WM> ' (B2) 

£C(u;,q) = Dio(w,q) = 0. (B4) 

(B5) 

The sign of is opposite the most common convention because we use a ( — h++) metric; 
but we have chosen the sign of I1t to correspond to common usage. The advanced propa- 
gator is given by complex conjugating these expressions. In what follows we drop the [ret] 
superscript from the self-energies. 
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It will be convenient for our purposes to define nonstandard normalization self-energies 

J2 



n T , = -n T 



n T = 



q 



q 2 — uS 1 



IT 



(B6) 



In terms of these, evaluating Eq. (Bl) gives 

1 



I t + U \ 



|g 2 + n L | 2 



+ 



(Apk-s)(Ap'k'-s) 



+ 



— ^— - — (p-k) (p-k)yj (Apk-s) (Ap'k'-s) cos 

((Apk-s){Apk'-s) cos 2 + q 2 (2uo 2 - Apk - Ap'k' 



|g 2 + n T | 



(B7) 



The analogous result for t channel exchange is 

1 



" / 2 i 2 \ " 



f 2 



|g 2 + n L | 2 



(App'+t)(Akk'+t) 



— ~ W2 - - (p+p')(k+k')^(App'+t)(Akk'+t) cos0 
(g 2 + n£)(g 2 + n T ) 

^— ((App'+t)(Akk'+t) cos 2 + g 2 (2w 2 + App + Akk')) 
- IIt 2 



|g 2 + n T | 



(B8) 



The variables in these expressions, such as 0, s, t, p, etc are defined in Appendix A. In 
particular the definition of in Eq. (B7) is the angle between the q, p plane and the q, p' 
plane, whereas in Eq. (B8) it is the angle between the q, p and q, k planes. 

Complete expressions for n L and I1 T were obtained by Weldon [27], and the usual ex- 
pressions for the hard thermal loop limit of the self-energies were extracted as a particular 
limit. The vacuum parts are simple 10 : 



n 



T, vac Tlh, vac 



12tt 2 



f 2 

q 



-log 



\u 2 -q 2 \ 
fi 2 



+ - +i7iQ{uj 2 -q 2 ) 



(B9) 



The imaginary part exists for timelike momenta and represents the vacuum rate of decay 
into a fermion pair. The full self-energy is the sum of a vacuum and a thermal part. The 
thermal parts of the self-energies are 



IIl, th — 



9 2 N^ 



H(u,q), n T ,th = 



9 2 NfT 2 ( 1 q 2 



G(u,q)} , (BIO) 



10 The earier (and published!) version is missing the 5/3 term here, an error which arose from an error in 
Weldon's paper [27] (unless his dimensionally regularized renormalization point a in his Eq. (A2) is meant 

to be /le 5 / 6 ). I thank Tony Rebhan for pointing out this error, which is rather significant at the largest 
values of g 2 Nf considered here. 



29 



where H(u, q) and G(u, q) are, 11 using the shorthand 2oj + = uj+q and 2c<j_ = uj— q, 




(Bll) 



The hard thermal loop limit means taking uj <C T, g <C T, and extracting the nonvanishing 
part. The values of these integrals, in this limit, depend only on u/q and are 



However at general uj/T, q/T, only a few of the integrals can be performed analytically and 
some have to be done numerically However, all required integrals can be done quickly and 
to high precision by numerical quadratures integration, and their evaluation does not limit 
either the accuracy or speed of any subsequent integrations. 

Note that, for timelike momenta, the thermal contributions to the imaginary parts of the 
self-energies are of opposite sign as the vacuum parts, and represent Pauli blocking of pair 
production. The thermal contribution to the imaginary part at spacelike momenta, absent 
in the vacuum theory, represents Landau damping. In the HTL limit only the imaginary 
part arising from Q(uj + )Q(— uj_) contributes. All real parts are even in uj and all imaginary 
parts are odd in uj. 



11 Our notation does not quite agree with [27]; we separate out the Debye mass squared g 2 NfT 2 /3 from the 
definitions of G and H, which are then pure numbers and have particularly simple small u, q limits. 



H UTL (uj/q) = 1 " f ^ -t + i— 8K)6(-,_) , 
2q uj- 2q 

G H TL(^/g) = 1 • 



(B12) 
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